Chapter 3 Spatial Autocorrelation Analysis

Considering the spatial lags of variables, variables in neighboring geographical units could be calculated, following the first order queen contiguity strategy, which classifies neighbors of geographical units by the sharing of either point or line segment borders. Followed by that, Global Moran’s I and Local Moran’s I statistics for charge points in London could be conducted.

3.1 Run Global Moran’s I

3.1.1 Function to calculate the centroids

# library(spdep)
# First calculate the centroids of all Wards in London

coordsWfun <- function(londonpointsyear){
  points_wards <- pointsjoinedfun(londonpointsyear)
  coordsW <- points_wards %>% 
    st_centroid() %>%
    st_geometry()
}

# test
coordsW <- coordsWfun(londonpoints2020)

## `summarise()` ungrouping output (override with `.groups` argument)

## Warning in st_centroid.sf(.): st_centroid assumes attributes are constant over geometries of x

plot(coordsW,axes=TRUE)

3.1.2 Function to create a neighbours list

#create a neighbours list
ward_nbfun <- function(londonpointsyear){
  points_wards <- pointsjoinedfun(londonpointsyear)
  ward_nb <- points_wards %>% 
    poly2nb(., queen=T)
  return(ward_nb)
}

# test
ward_nb20 <- ward_nbfun(londonpoints2020)

## `summarise()` ungrouping output (override with `.groups` argument)

plot(ward_nb20, st_geometry(coordsW), col="red")
# add a map underneath
points_wards <- pointsjoinedfun(londonpoints2020) 

## `summarise()` ungrouping output (override with `.groups` argument)

plot(points_wards$geometry, add=T)

3.1.3 Function to create spatial weights

#create a spatial weights object from these weights
ward_lwfun <- function(londonpointsyear){
  ward_lw <- ward_nbfun(londonpointsyear) %>% 
    nb2listw(., style="C")
  return(ward_lw)
}

# test
ward_lw20 <- ward_lwfun(londonpoints2020)

## `summarise()` ungrouping output (override with `.groups` argument)

head(ward_lw20$neighbours)

## [[1]]
## [1]   6   7  10 462 468 482
## 
## [[2]]
## [1]  5  8  9 11 12 13 16
## 
## [[3]]
## [1]  10  11  12  15 480 483
## 
## [[4]]
## [1]  17 281 291 470 473 481
## 
## [[5]]
## [1]   2   9  16 281 283 290
## 
## [[6]]
## [1]  1  7  8 10 11 14

3.1.4 Function to run Global Moran’s I

I_ward_globalfun <- function(londonpointsyear){
  ward_lw <- ward_lwfun(londonpointsyear)
  
  I_ward_global <- londonpointsyear %>% 
    pointsjoinedfun(.) %>% 
    pull(density) %>%
    as.vector()%>% # Converts a distributed matrix into a non-distributed vector
    moran.test(., ward_lw) # density <-> spatial (similar or not)
  return(I_ward_global)
}

3.1.5 Global Moran’s I for 2018, 2019, 2020

I_ward_global18 <- I_ward_globalfun(londonpoints2018)

## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)

I_ward_global19 <- I_ward_globalfun(londonpoints2019)

## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)

I_ward_global20 <- I_ward_globalfun(londonpoints2020)

## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)

I_ward_global18

## 
##  Moran I test under randomisation
## 
## data:  .  
## weights: ward_lw    
## 
## Moran I statistic standard deviate = 11.54, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##      0.2589508028     -0.0016025641      0.0005097484

I_ward_global19

## 
##  Moran I test under randomisation
## 
## data:  .  
## weights: ward_lw    
## 
## Moran I statistic standard deviate = 23.078, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##      0.5238081425     -0.0016025641      0.0005183155

I_ward_global20

## 
##  Moran I test under randomisation
## 
## data:  .  
## weights: ward_lw    
## 
## Moran I statistic standard deviate = 25.345, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##       0.579937112      -0.001602564       0.000526483

3.2 Run Local Moran’s I

3.2.1 Function to run Local Moran’s I

ward_lw20 <- ward_lwfun(londonpoints2020)

## `summarise()` ungrouping output (override with `.groups` argument)

  I_ward_local_count <- londonpoints2020 %>% 
    pointsjoinedfun(.) %>% 
    pull(chargepointcount) %>%
    as.vector()%>% 
    localmoran(., ward_lw20) %>% 
    as_tibble()

## `summarise()` ungrouping output (override with `.groups` argument)

  I_ward_local_density <- londonpoints2020 %>% 
    pointsjoinedfun(.) %>% 
    pull(density) %>%
    as.vector()%>% 
    localmoran(., ward_lw20) %>% 
    as_tibble()

## `summarise()` ungrouping output (override with `.groups` argument)

slice_head(I_ward_local_count, n=5)

## # A tibble: 5 x 5
##   Ii         E.Ii         Var.Ii     Z.Ii        `Pr(z > 0)`
##   <localmrn> <localmrn>   <localmrn> <localmrn>  <localmrn> 
## 1 -0.0152550 -0.001633048 0.1694194  -0.03309466 0.5132004  
## 2 -0.1351878 -0.001905222 0.1973404  -0.30003036 0.6179230  
## 3  0.4127932 -0.001633048 0.1694194   1.00685228 0.1570029  
## 4  0.3218993 -0.001633048 0.1694194   0.78602481 0.2159265  
## 5 -0.1455548 -0.001633048 0.1694194  -0.34965929 0.6367028

slice_head(I_ward_local_density, n=5)

## # A tibble: 5 x 5
##   Ii          E.Ii         Var.Ii     Z.Ii        `Pr(z > 0)`
##   <localmrn>  <localmrn>   <localmrn> <localmrn>  <localmrn> 
## 1 -0.02311287 -0.001633048 0.1675912  -0.05246928 0.5209226  
## 2  0.09282148 -0.001905222 0.1952145   0.21439584 0.4151192  
## 3  0.25021517 -0.001633048 0.1675912   0.61519570 0.2692127  
## 4  0.27698598 -0.001633048 0.1675912   0.68058939 0.2480657  
## 5  0.07334936 -0.001633048 0.1675912   0.18316134 0.4273357

# use the local moran function to generate I for each ward in the city

I_ward_localfun <- function(londonpointsyear){
  ward_lw <- ward_lwfun(londonpointsyear)
  
  I_ward_local_count <- londonpointsyear %>% 
    pointsjoinedfun(.) %>% 
    pull(chargepointcount) %>%
    as.vector()%>% 
    localmoran(., ward_lw) %>% 
    as_tibble()
  I_ward_local_density <- londonpointsyear %>% 
    pointsjoinedfun(.) %>% 
    pull(density) %>%
    as.vector()%>% 
    localmoran(., ward_lw) %>% 
    as_tibble()
  
  tm_local_moran <- londonpointsyear %>% 
    pointsjoinedfun(.) %>%
    mutate(charge_count_I = as.numeric(I_ward_local_count$Ii))%>%
    mutate(charge_count_Iz =as.numeric(I_ward_local_count$Z.Ii))%>%
    mutate(density_I =as.numeric(I_ward_local_density$Ii))%>%
    mutate(density_Iz =as.numeric(I_ward_local_density$Z.Ii))
  
  return(tm_local_moran)
}

3.2.2 Function to run tmap for Local Moran’s I

tmap_local_moranfun <- function(londonpointsyear, number, mode){
  breaks1<-c(-1000,-2.58,-1.96,-1.65,1.65,1.96,2.58,1000)
  MoranColours<- rev(brewer.pal(8, "RdBu"))
  tm_local_moran <- I_ward_localfun(londonpointsyear)
  tmap_mode(mode)
  tmap <- tm_shape(tm_local_moran) +
    tm_polygons("charge_count_Iz",
        style="fixed",
        breaks=breaks1,
        palette=MoranColours,
        midpoint=NA)+
    tm_legend(show=FALSE)+
    tm_layout(frame=FALSE)+
    tm_credits(number, position=c(0,0.85), size=1)
  return(tmap)
}

tmap_moranlegendfun <- function(londonpointsyear){
  breaks1<-c(-1000,-2.58,-1.96,-1.65,1.65,1.96,2.58,1000)
  MoranColours<- rev(brewer.pal(8, "RdBu"))
  tm_local_moran <- I_ward_localfun(londonpointsyear)
  legend <- tm_shape(tm_local_moran) +
    tm_polygons("charge_count_Iz",
                breaks=breaks1,
                palette=MoranColours,
                title="Local Moran's I, Charge Points in London") +
    tm_scale_bar(position=c(0.1,0.04), text.size=0.6)+
    tm_compass(north=0, position=c(0.65,0.2))+
   
    tm_layout(title = "Charge Points Local Moran's I", 
              legend.title.size=1,
              legend.text.size = 0.6,
              legend.only = TRUE, 
              legend.position=c(0.1,0.25),asp=0.1)
  
  return(legend)
}

tm_local_moran18 <- tmap_local_moranfun(londonpoints2018, "a)", "plot")

## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)

## tmap mode set to plotting

tm_local_moran19 <- tmap_local_moranfun(londonpoints2019, "b)", "plot")

## `summarise()` ungrouping output (override with `.groups` argument)

## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)

## tmap mode set to plotting

tm_local_moran20 <- tmap_local_moranfun(londonpoints2020, "c)", "plot")

## `summarise()` ungrouping output (override with `.groups` argument)

## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)

## tmap mode set to plotting

moran_legend <- tmap_moranlegendfun(londonpoints2020)

## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)

local_moran <- tmap_arrange(tm_local_moran18, 
                            tm_local_moran19, 
                            tm_local_moran20, 
                            moran_legend, ncol=2)
local_moran

## Variable(s) "charge_count_Iz" contains positive and negative values, so midpoint is set to 0. Set midpoint = NA to show the full spectrum of the color palette.

tmap_save(local_moran, "pic/Local Moran's I, Charge Points in London.png",width=7, height=4)

## Variable(s) "charge_count_Iz" contains positive and negative values, so midpoint is set to 0. Set midpoint = NA to show the full spectrum of the color palette.

## Map saved to /Users/apple/OneDrive - University College London/Module assessment/CASA_project/casa_project/pic/Local Moran's I, Charge Points in London.png

## Resolution: 2100 by 1200 pixels

## Size: 7 by 4 inches (300 dpi)

tm_local_moran18v <- tmap_local_moranfun(londonpoints2018, "a)", "view")

## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)

## tmap mode set to interactive viewing

tm_local_moran19v <- tmap_local_moranfun(londonpoints2019, "b)", "view")

## `summarise()` ungrouping output (override with `.groups` argument)

## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)

## tmap mode set to interactive viewing

tm_local_moran20v <- tmap_local_moranfun(londonpoints2020, "c)", "view")

## `summarise()` ungrouping output (override with `.groups` argument)

## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)

## tmap mode set to interactive viewing

moran_legend <- tmap_moranlegendfun(londonpoints2020)

## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)

local_moran <- tmap_arrange(tm_local_moran18v, 
                            tm_local_moran19v, 
                            tm_local_moran20v, ncol=1)
local_moran

## Credits not supported in view mode.

## Credits not supported in view mode.
## Credits not supported in view mode.

3.2.3 Function to run Gi*

Gi_ward_local_densityfun <- function(londonpointsyear){
  ward_lw <- ward_lwfun(londonpointsyear)
  
  Gi_ward_local_density <- londonpointsyear %>% 
    pointsjoinedfun(.) %>% 
    pull(density) %>%
    as.vector()%>%
    localG(., ward_lw)
  
  density_Gi <- londonpointsyear %>% 
    pointsjoinedfun(.) %>% 
    mutate(density_G = as.numeric(Gi_ward_local_density))
  
  return(density_Gi)
}

3.2.4 Function to run tmap for Gi*

tmap_Gifun <- function(londonpointsyear, number, mode){
  breaks1<-c(-1000,-2.58,-1.96,-1.65,1.65,1.96,2.58,1000)
  GIColours<- rev(brewer.pal(8, "RdBu"))
  density_Gi <- Gi_ward_local_densityfun(londonpointsyear)
  tmap_mode(mode)
  tmap <- tm_shape(density_Gi) +
    tm_polygons("density_G",
                breaks=breaks1,
                palette=GIColours,
                title="Gi*, Charge Points in London") +
    tm_legend(show=FALSE)+
    tm_layout(frame=FALSE)+
    tm_credits(number, position=c(0,0.85), size=1)
  
  return(tmap)
}

G1 <- tmap_Gifun(londonpoints2018, "a)", "plot")

## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)

## tmap mode set to plotting

G2 <- tmap_Gifun(londonpoints2019, "b)", "plot")

## `summarise()` ungrouping output (override with `.groups` argument)

## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)

## tmap mode set to plotting

G3 <- tmap_Gifun(londonpoints2020, "c)", "plot")

## `summarise()` ungrouping output (override with `.groups` argument)

## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)

## tmap mode set to plotting

tmap_Gilegendfun <- function(londonpointsyear){
  breaks1<-c(-1000,-2.58,-1.96,-1.65,1.65,1.96,2.58,1000)
  GIColours<- rev(brewer.pal(8, "RdBu"))
  density_Gi <- Gi_ward_local_densityfun(londonpointsyear)
  legend <- tm_shape(density_Gi) +
    tm_polygons("density_G",
                breaks=breaks1,
                palette=GIColours,
                title="Gi*, Charge Points in London") +
    tm_scale_bar(position=c(0.1,0.04), text.size=0.6)+
    tm_compass(north=0, position=c(0.65,0.2))+
   
    tm_layout(title = "Charge Points Gi* Density", 
              legend.title.size=1,
              legend.text.size = 0.6,
              legend.only = TRUE, 
              legend.position=c(0.1,0.25),asp=0.1)
  
  return(legend)
}

Gi_legend <- tmap_Gilegendfun(londonpoints2020)

## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)
## `summarise()` ungrouping output (override with `.groups` argument)

Gi <- tmap_arrange(G1, G2, G3, Gi_legend, ncol=2)
Gi

## Variable(s) "density_G" contains positive and negative values, so midpoint is set to 0. Set midpoint = NA to show the full spectrum of the color palette.

## Variable(s) "density_G" contains positive and negative values, so midpoint is set to 0. Set midpoint = NA to show the full spectrum of the color palette.
## Variable(s) "density_G" contains positive and negative values, so midpoint is set to 0. Set midpoint = NA to show the full spectrum of the color palette.
## Variable(s) "density_G" contains positive and negative values, so midpoint is set to 0. Set midpoint = NA to show the full spectrum of the color palette.

tmap_save(local_moran, "pic/Gi, Charge Points in London.png",width=7, height=4)

## Map saved to /Users/apple/OneDrive - University College London/Module assessment/CASA_project/casa_project/pic/Gi, Charge Points in London.png

## Resolution: 2100 by 1200 pixels

## Size: 7 by 4 inches (300 dpi)